Introduction

Raw data

The raw data is shown below. There were 196 rows and 42 columns, consisting of 28 patients taken over 7 timepoints. However, there is a significant amount of missing data, resulting in only 120 usable datapoints.

missing.ix <- Reduce(intersect, apply(df[,14:42], 2, function(x) which(is.na(x))))
df.raw <- df[-missing.ix,]
df.raw <- df.raw[order(df.raw$PatientID),]
rownames(df.raw) <- 1:nrow(df.raw)
kable(df.raw, 
      digits = 3,
      row.names = T,
      caption = "Raw Data"
) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12) %>%
    scroll_box(width = "100%", height = "300px") 
Raw Data
PatientID Age AgeGreater60 Sex LowIntermacs InterMACS RVAD Sensitized VAD Indication Device Type Survival Outcome Time num Total PBMC num lymph lymph live lymph CD3 of live lymph CD19 of live lymph CD19+CD27- CD19+CD27+ CD27+38++plasma blasts CD27-38++ transitional CD27-IgD+ mature naive CD27+IgD- switched memory CD27-IgD- switched memory CD27+IgD+ unswitched memory CD27+IgD-IgM+ switched memory CD27+IgD+IgM+ nonswitched memory CD19+27+IgG+IgM- memory CD19+24dim38dim naive mature CD19+24+38++transitional CD19CD24hiCD38-memory CD19+27-38+CD5+transitionals CD19+CD268+ CD268 of +27-38++transitional CD19+CD11b+ CD19+CD5+ CD19+CD27+CD24hi CD19+CD5+CD24hi CD19+CD5+CD11b+ CD19+27+IgD-38++IgG ASC
1 1 65 older Male High 2 No NA BTT HMII alive Alive s/p OHT 0 169154 35496 20.98 99.53 25.37 19.82 86.02 13.98 1.26 2.54 57.81 11.74 28.06 2.39 21.74 14.36 0.72 81.05 2.86 13.60 0.58 95.22 84.83 10.47 4.33 8.50 1.39 3.13 1.09
2 1 65 older Male High 2 No NA BTT HMII alive Alive s/p OHT 1 63082 9915 15.72 99.26 32.45 26.26 90.06 9.94 0.89 2.09 54.27 7.70 35.71 2.32 31.10 20.47 1.57 84.06 2.59 10.79 0.26 97.02 92.59 7.50 3.13 6.42 1.32 2.44 1.49
3 1 65 older Male High 2 No NA BTT HMII alive Alive s/p OHT 3 75921 21721 28.61 99.31 24.86 33.01 86.38 13.62 1.54 1.26 45.74 10.98 40.44 2.84 31.80 18.04 1.15 86.21 4.82 6.56 0.24 90.41 73.33 10.00 7.25 10.12 3.22 5.67 2.05
4 1 65 older Male High 2 No NA BTT HMII alive Alive s/p OHT 8 213808 36002 16.84 99.19 35.95 14.55 67.94 32.06 4.06 2.17 27.92 28.92 39.64 3.52 20.31 9.10 1.33 68.69 7.10 20.42 1.88 87.38 44.25 10.78 11.12 19.15 5.54 7.35 1.07
5 1 65 older Male High 2 No NA BTT HMII alive Alive s/p OHT 21 106970 31867 29.79 99.18 42.36 15.02 84.69 15.31 1.94 2.11 55.41 13.14 29.13 2.32 20.22 12.88 2.08 83.21 5.20 9.50 1.10 94.44 73.00 9.04 6.99 9.22 2.80 4.78 1.44
6 2 81 older Male Low 3 No NA DT HMII dead Died 0 119377 77509 64.93 99.61 58.50 7.99 62.67 37.33 9.10 4.82 36.73 31.90 25.62 5.74 11.47 8.95 0.40 69.54 7.69 9.41 2.61 76.44 57.91 23.29 20.76 15.39 5.42 14.90 2.60
7 2 81 older Male Low 3 No NA DT HMII dead Died 3 96940 57658 59.48 99.70 55.52 9.21 72.15 27.85 3.85 3.82 48.41 21.75 23.49 6.35 12.10 16.04 0.28 75.30 3.72 6.65 2.11 80.52 67.33 15.80 15.23 6.71 1.32 8.86 0.26
8 2 81 older Male Low 3 No NA DT HMII dead Died 5 154373 73969 47.92 99.65 45.18 8.77 67.36 32.64 6.59 2.77 39.10 26.28 27.95 6.67 17.89 15.71 0.39 71.57 6.26 9.88 1.78 81.86 72.07 20.12 20.17 14.80 5.99 14.65 3.71
9 2 81 older Male Low 3 No NA DT HMII dead Died 8 155610 63375 40.73 99.51 38.42 12.40 71.69 28.31 4.36 2.62 37.82 23.43 33.50 5.24 17.37 14.92 0.60 65.04 2.80 11.23 2.54 83.55 37.56 15.55 15.45 9.55 2.66 10.91 0.28
10 2 81 older Male Low 3 No NA DT HMII dead Died 14 244245 115109 47.13 99.47 47.40 9.37 67.50 32.50 3.99 2.08 50.08 26.01 17.22 6.68 15.34 15.25 0.23 77.86 3.21 10.49 1.48 82.13 59.19 21.45 19.77 14.45 2.57 14.94 0.94
11 3 58 younger Male Low 3 Yes No BTT HMII alive Alive s/p OHT 1 246811 86529 35.06 93.49 67.96 17.48 81.66 18.34 2.05 6.85 51.33 13.69 29.91 5.08 15.68 5.94 0.19 61.95 0.76 13.51 4.07 92.62 34.06 6.90 12.07 8.00 2.54 3.85 0.26
12 3 58 younger Male Low 3 Yes No BTT HMII alive Alive s/p OHT 3 406771 90343 22.21 97.34 49.28 17.02 75.04 24.96 4.47 18.36 28.45 21.86 46.07 3.62 24.97 4.33 0.13 47.41 1.64 19.70 18.15 72.76 8.08 25.69 28.35 13.95 4.66 19.20 0.41
13 3 58 younger Male Low 3 Yes No BTT HMII alive Alive s/p OHT 21 330191 81636 24.72 94.57 45.14 14.70 70.51 29.49 8.94 16.23 38.92 23.53 30.83 6.72 23.93 7.07 0.18 46.63 2.35 16.12 8.59 74.45 39.90 24.19 21.18 12.62 3.26 12.93 0.38
14 4 67 older Male Low 3 No Yes BTT HVAD alive Alive s/p OHT 0 309890 41896 13.52 97.03 47.13 8.08 77.31 22.69 2.74 3.29 59.03 19.80 19.37 1.80 21.10 4.46 19.79 71.98 2.41 17.39 0.20 95.80 75.00 5.79 1.89 16.45 0.91 0.67 0.77
15 4 67 older Male Low 3 No Yes BTT HVAD alive Alive s/p OHT 1 274228 15653 5.71 97.85 35.76 28.87 51.45 48.55 2.60 2.13 30.85 43.13 22.77 3.26 20.09 3.92 30.51 57.58 1.31 28.20 0.18 91.23 59.57 12.08 2.04 34.44 0.38 0.16 0.16
16 4 67 older Male Low 3 No Yes BTT HVAD alive Alive s/p OHT 3 591939 57589 9.73 96.76 38.20 42.54 75.71 24.29 2.14 1.33 53.03 21.26 24.16 1.55 25.50 4.02 20.27 67.18 1.14 21.71 0.11 95.28 66.98 3.79 2.43 17.70 0.65 0.84 0.18
17 4 67 older Male Low 3 No Yes BTT HVAD alive Alive s/p OHT 5 189634 46711 24.63 97.74 34.84 37.32 75.30 24.70 0.70 0.55 46.97 22.05 30.28 0.70 23.67 1.56 31.44 73.33 0.75 18.59 0.00 84.61 86.17 3.44 2.11 18.68 0.94 0.63 0.69
18 4 67 older Male Low 3 No Yes BTT HVAD alive Alive s/p OHT 14 182319 40475 22.20 96.28 43.36 11.77 69.34 30.66 1.59 0.89 52.64 27.32 18.56 1.48 22.79 2.88 28.61 75.47 0.33 14.54 0.06 95.90 53.66 7.81 3.27 21.04 1.57 1.50 0.72
19 5 65 older Female High 2 No Yes BTT HMII alive Alive s/p OHT 0 1073919 336538 31.34 100.00 63.29 14.80 93.05 6.95 0.73 8.05 75.49 4.55 19.16 0.80 7.25 8.78 28.04 88.72 2.51 5.18 1.02 90.49 68.18 2.32 4.09 1.69 3.48 1.27 5.78
20 5 65 older Female High 2 No Yes BTT HMII alive Alive s/p OHT 1 270301 97229 35.97 91.63 18.43 61.67 97.68 2.32 0.20 5.85 69.45 1.63 28.19 0.73 27.57 21.43 17.05 84.32 7.26 2.73 0.84 99.44 97.60 0.54 2.64 1.51 3.63 0.29 3.40
21 5 65 older Female High 2 No Yes BTT HMII alive Alive s/p OHT 3 499016 78019 15.63 89.41 25.83 51.38 92.98 7.02 0.32 4.66 51.34 6.31 41.59 0.76 36.49 5.96 22.81 73.48 6.73 7.93 0.33 98.78 96.52 1.31 2.06 5.19 2.68 0.59 2.45
22 5 65 older Female High 2 No Yes BTT HMII alive Alive s/p OHT 8 115825 32513 28.07 91.04 44.61 15.76 80.52 19.48 4.16 7.82 49.91 18.37 30.52 1.20 14.87 3.33 30.19 61.08 8.70 12.56 1.58 89.80 76.71 2.87 6.07 11.12 8.08 1.59 1.17
23 5 65 older Female High 2 No Yes BTT HMII alive Alive s/p OHT 21 407080 97651 23.99 91.29 66.96 12.20 84.66 15.34 2.88 5.60 61.89 13.73 22.65 1.74 20.16 7.23 29.39 68.49 5.28 10.60 1.17 92.50 54.35 3.48 6.53 9.61 8.34 1.85 2.83
24 6 43 younger Male High 2 Yes No BTT HMII alive Alive s/p OHT 0 515760 266115 51.60 96.30 69.66 8.57 72.22 27.78 6.26 4.77 43.09 26.01 28.26 2.63 16.67 2.67 28.25 50.35 0.41 12.78 4.47 84.76 19.47 9.56 4.55 8.29 0.87 2.18 0.75
25 6 43 younger Male High 2 Yes No BTT HMII alive Alive s/p OHT 8 539521 223872 41.49 96.72 60.02 12.63 80.13 19.87 1.66 2.48 43.56 19.52 35.24 1.68 20.46 3.67 20.54 55.76 0.18 10.80 2.43 87.82 13.25 6.55 5.96 4.65 0.95 2.51 0.65
26 6 43 younger Male High 2 Yes No BTT HMII alive Alive s/p OHT 21 564317 216312 38.33 97.40 55.80 4.83 74.19 25.81 5.35 5.77 37.90 25.72 34.76 1.62 17.99 2.18 31.75 61.61 0.77 13.44 6.97 76.55 6.97 10.34 11.14 4.46 0.94 5.14 2.43
27 7 49 younger Male Low 3 No No BTT HVAD alive Alive s/p OHT 0 89251 65181 73.03 76.76 84.53 1.15 74.65 25.35 1.91 1.39 55.03 22.40 19.62 2.95 15.65 6.80 25.85 70.49 0.87 23.26 0.23 83.16 25.00 8.85 5.03 11.63 1.56 2.95 3.08
28 7 49 younger Male Low 3 No No BTT HVAD alive Alive s/p OHT 3 375167 234817 62.59 98.46 60.29 10.18 87.72 12.28 0.46 0.48 68.54 9.90 19.03 2.54 30.58 19.69 15.54 72.77 0.13 20.35 0.03 92.56 36.28 3.75 2.70 9.17 1.19 1.27 0.63
29 7 49 younger Male Low 3 No No BTT HVAD alive Alive s/p OHT 5 256306 162176 63.27 93.82 68.43 8.45 85.45 14.55 0.70 0.53 71.99 11.98 13.37 2.66 32.41 17.96 12.46 69.74 0.41 24.63 0.02 95.72 44.12 4.47 2.23 12.99 1.30 1.20 0.77
30 7 49 younger Male Low 3 No No BTT HVAD alive Alive s/p OHT 8 463136 330714 71.41 98.79 70.51 7.30 87.51 12.49 1.78 0.58 69.60 10.61 17.82 1.97 29.83 14.55 12.59 75.03 0.15 21.61 0.07 93.59 30.94 4.05 3.07 10.00 1.08 1.35 0.62
31 8 43 younger Male Low 3 No No BTT HMII alive Alive s/p OHT 0 92221 27988 30.35 85.45 55.78 1.51 85.56 14.44 0.28 0.28 55.28 4.72 39.72 0.28 0.00 1.89 1.89 89.44 0.00 2.78 0.00 15.00 0.00 3.33 0.83 0.28 0.28 3.33 0.00
32 8 43 younger Male Low 3 No No BTT HMII alive Alive s/p OHT 1 273746 46755 17.08 93.25 50.14 7.57 88.27 11.73 0.09 2.12 54.52 3.36 42.00 0.12 8.40 6.87 0.00 92.00 0.12 3.45 0.00 3.94 11.43 1.12 0.64 2.76 0.33 2.15 0.00
33 8 43 younger Male Low 3 No No BTT HMII alive Alive s/p OHT 3 156053 71752 45.98 78.38 35.43 27.12 78.57 21.43 0.24 0.60 58.56 5.80 18.33 17.31 12.67 72.96 0.53 94.06 0.49 5.11 0.28 22.79 27.47 0.81 5.53 0.01 0.29 0.19 0.00
34 8 43 younger Male Low 3 No No BTT HMII alive Alive s/p OHT 8 156053 71752 45.98 78.38 35.43 27.12 78.57 21.43 0.24 0.60 58.56 5.80 18.33 17.31 12.67 72.96 0.53 94.06 0.49 5.11 0.28 22.79 27.47 0.81 5.53 0.01 0.29 0.19 0.00
35 9 60 younger Female High 2 Yes Yes BTT PVAD dead Died 3 116462 41076 35.27 95.44 33.79 25.24 84.51 15.49 2.92 0.99 45.34 13.19 39.12 2.36 17.56 8.75 19.73 79.35 5.17 10.80 0.71 90.04 29.59 10.16 6.81 8.21 3.37 4.86 4.47
36 9 60 younger Female High 2 Yes Yes BTT PVAD dead Died 5 222650 85649 38.47 98.40 32.09 19.82 77.54 22.46 5.02 1.10 43.96 19.20 33.50 3.34 17.47 9.77 15.63 71.89 6.68 14.08 0.70 86.49 27.87 13.04 9.84 12.30 4.35 7.98 6.39
37 9 60 younger Female High 2 Yes Yes BTT PVAD dead Died 8 102274 63057 61.65 97.73 27.89 20.66 73.75 26.25 7.33 3.45 35.79 21.92 37.73 4.56 9.97 8.20 16.96 63.96 2.40 11.16 3.07 73.97 5.92 16.69 17.77 11.74 4.90 12.19 6.09
38 9 60 younger Female High 2 Yes Yes BTT PVAD dead Died 14 133530 81665 61.16 96.46 22.03 16.15 35.30 64.70 13.74 2.33 8.48 52.28 26.37 12.87 5.29 5.79 19.15 47.48 5.65 17.83 5.43 39.88 6.08 52.08 62.11 37.45 21.95 50.15 3.84
39 10 36 younger Male High 2 Yes NA BTT PVAD alive Alive s/p OHT 0 481928 162052 33.63 99.70 64.73 5.05 3.91 96.09 8.89 1.34 0.72 94.99 3.14 1.15 4.07 0.84 0.32 90.28 2.84 1.58 1.27 6.82 9.17 14.61 90.11 4.95 5.71 13.01 0.51
40 10 36 younger Male High 2 Yes NA BTT PVAD alive Alive s/p OHT 8 38113 21486 56.37 97.86 58.44 7.88 20.47 79.53 8.09 1.69 1.21 75.00 18.96 4.83 10.07 2.59 0.53 82.19 2.29 9.00 3.59 15.16 0.00 23.07 63.95 4.77 6.52 16.91 0.40
41 10 36 younger Male High 2 Yes NA BTT PVAD alive Alive s/p OHT 21 80991 19210 23.72 99.00 50.61 3.50 12.33 87.67 14.74 4.81 1.05 85.41 11.13 2.41 8.48 0.52 1.56 80.00 8.27 3.76 7.41 13.68 0.00 30.83 77.74 10.83 14.14 19.85 2.81
42 11 61 older Male High 2 No No BTT HMII alive Alive s/p OHT 0 411616 119210 28.96 90.82 73.63 5.88 74.30 25.70 2.80 6.46 57.94 19.20 18.72 4.13 19.10 15.42 15.30 53.51 13.78 25.19 0.21 47.23 53.77 5.75 5.92 22.63 3.46 2.55 0.74
43 11 61 older Male High 2 No No BTT HMII alive Alive s/p OHT 1 430509 79492 18.46 89.31 52.10 20.45 87.13 12.87 1.23 4.76 58.11 9.15 30.55 2.20 29.86 18.05 10.37 63.63 10.23 17.43 0.06 74.05 72.50 2.36 2.98 10.68 1.72 1.03 2.89
44 11 61 older Male High 2 No No BTT HMII alive Alive s/p OHT 3 395713 119670 30.24 97.90 23.79 6.27 82.17 17.83 0.79 1.99 49.95 13.01 33.91 3.13 28.36 16.83 5.00 56.33 7.13 26.73 0.02 79.50 82.88 1.96 5.55 16.30 3.61 1.88 0.42
45 11 61 older Male High 2 No No BTT HMII alive Alive s/p OHT 8 329840 73468 22.27 98.63 70.61 3.73 51.50 48.50 8.32 3.77 33.46 41.29 20.92 4.33 25.19 5.89 24.07 46.40 3.73 34.64 0.07 64.14 28.43 10.35 9.61 38.78 4.62 5.55 2.91
46 11 61 older Male High 2 No No BTT HMII alive Alive s/p OHT 14 350294 88823 25.36 97.63 74.37 4.37 73.35 26.65 4.91 4.83 56.54 21.32 19.84 2.30 34.68 11.25 18.37 61.98 5.78 20.55 0.04 63.91 47.54 6.81 6.17 21.77 2.90 2.74 3.00
47 12 38 younger Male High 2 Yes Yes BTT TAH dead Died s/p OHT 0 128451 60927 47.43 96.82 60.77 4.70 77.26 22.74 2.17 1.95 55.58 18.40 21.54 4.47 27.86 18.89 15.33 69.90 0.94 24.40 0.14 70.73 1.85 8.01 8.66 20.21 5.45 3.32 0.39
48 12 38 younger Male High 2 Yes Yes BTT TAH dead Died s/p OHT 1 142796 22573 15.81 80.67 47.31 2.28 63.37 36.63 9.16 6.51 26.02 33.49 37.11 3.37 38.13 7.50 13.13 53.25 5.78 26.27 1.57 40.24 0.00 27.71 22.17 29.16 18.80 15.66 8.11
49 12 38 younger Male High 2 Yes Yes BTT TAH dead Died s/p OHT 3 157167 38008 24.18 93.47 34.86 3.95 69.71 30.29 3.92 3.92 33.64 26.30 35.78 4.28 39.54 11.49 14.94 66.86 2.85 24.31 1.04 64.50 7.27 21.60 20.53 27.66 14.54 14.18 2.12
50 12 38 younger Male High 2 Yes Yes BTT TAH dead Died s/p OHT 5 121823 23097 18.96 93.48 2.51 2.87 99.52 0.48 0.00 0.00 0.97 0.48 98.55 0.00 50.00 25.00 25.00 94.35 0.00 2.58 0.00 0.16 0.00 17.45 87.40 0.00 0.16 18.26 0.00
51 12 38 younger Male High 2 Yes Yes BTT TAH dead Died s/p OHT 8 196515 88579 45.07 85.54 17.17 3.28 59.38 40.62 3.21 1.41 35.44 35.28 23.87 5.42 62.85 12.32 3.71 66.05 4.38 23.22 0.27 74.53 31.43 25.51 28.12 39.78 25.79 21.57 0.79
52 13 61 older Male High 1 Yes No BTT TAH alive Alive s/p OHT 1 87827 15751 17.93 76.38 31.60 12.17 77.46 22.54 6.35 2.53 50.34 20.22 27.05 2.39 19.68 7.62 19.05 75.55 7.38 9.49 0.71 94.26 89.19 11.68 17.96 18.78 6.97 10.86 12.46
53 13 61 older Male High 1 Yes No BTT TAH alive Alive s/p OHT 3 55539 13608 24.50 79.00 28.02 15.78 81.90 18.10 6.13 2.42 66.92 14.33 14.98 3.77 20.81 17.45 18.12 78.01 7.08 8.61 0.72 94.04 70.73 11.32 16.69 15.21 5.78 10.50 16.39
54 13 61 older Male High 1 Yes No BTT TAH alive Alive s/p OHT 8 295005 36785 12.47 94.58 36.55 9.14 86.32 13.68 2.45 1.82 37.01 12.52 49.18 1.29 16.24 4.71 12.47 8.36 1.73 6.64 0.81 84.78 46.55 12.52 9.28 5.94 1.73 4.47 2.49
55 13 61 older Male High 1 Yes No BTT TAH alive Alive s/p OHT 14 60695 20471 33.73 94.91 6.93 7.42 83.97 16.03 3.40 0.69 24.36 14.71 59.33 1.60 18.43 5.07 11.52 54.82 2.91 16.10 0.50 82.44 30.00 51.08 13.25 10.41 3.68 11.31 6.06
56 14 39 younger Female Low 3 Yes No BTT HMII dead Died post OHT 0 52016 30920 59.44 95.77 55.57 24.62 58.73 41.27 2.14 0.97 53.01 5.05 5.45 36.50 9.39 81.29 2.88 81.48 5.65 9.02 0.48 97.75 71.83 0.67 22.89 0.34 1.78 0.49 NA
57 14 39 younger Female Low 3 Yes No BTT HMII dead Died post OHT 1 75516 48801 64.62 92.13 52.90 18.64 70.39 29.61 2.92 2.21 62.18 5.13 7.98 24.70 13.79 69.54 3.19 84.14 3.64 6.36 0.41 94.93 43.78 0.57 15.61 0.23 0.91 0.36 NA
58 14 39 younger Female Low 3 Yes No BTT HMII dead Died post OHT 3 250612 130454 52.05 95.75 60.05 3.29 46.04 53.96 1.34 3.89 33.08 32.47 8.29 26.15 3.56 34.61 0.09 76.71 0.39 13.20 0.28 39.79 14.38 2.26 1.56 10.45 2.89 3.04 0.08
59 14 39 younger Female Low 3 Yes No BTT HMII dead Died post OHT 8 212478 120477 56.70 85.12 48.36 12.51 79.57 20.43 1.35 2.14 58.73 12.87 18.44 9.96 5.63 34.32 0.11 86.74 0.41 6.88 0.15 35.30 53.09 1.23 0.75 2.36 1.04 1.47 NA
60 15 70 older Female High 2 No NA DT HMII alive Alive 0 86135 70306 81.62 93.71 74.50 1.59 61.74 38.26 1.91 2.67 33.97 14.50 26.81 24.71 20.67 62.47 9.03 73.00 2.00 24.71 1.27 20.23 71.43 0.48 12.79 0.38 3.24 0.48 4.60
61 15 70 older Female High 2 No NA DT HMII alive Alive 1 81004 41334 51.03 94.06 53.64 2.99 61.70 38.30 3.10 2.75 39.07 13.34 22.12 25.47 12.58 63.11 15.99 78.92 2.93 17.38 3.35 23.24 71.88 0.52 11.70 0.09 1.29 0.52 9.83
62 15 70 older Female High 2 No NA DT HMII alive Alive 5 90024 39434 43.80 95.32 68.83 5.62 71.59 28.41 0.85 1.66 35.23 13.02 35.80 15.96 22.94 52.37 5.85 76.42 1.14 21.83 0.62 14.54 48.57 0.38 4.69 0.00 0.80 0.38 0.62
63 15 70 older Female High 2 No NA DT HMII alive Alive 14 106791 80532 75.41 88.82 79.92 1.57 62.32 37.68 1.70 4.02 28.04 22.86 34.11 15.00 23.92 39.41 15.03 61.61 3.66 34.55 1.97 20.36 73.33 0.80 8.57 0.27 1.79 0.98 0.71
64 16 63 older Male High 1 No NA BTT HMII dead Died 3 366713 77282 21.07 98.61 66.55 12.44 93.65 6.35 0.53 0.60 62.22 4.68 32.38 0.72 19.94 9.89 28.23 83.18 0.51 14.38 0.30 94.72 54.39 2.07 11.33 2.86 3.40 1.48 3.58
65 16 63 older Male High 1 No NA BTT HMII dead Died 5 262251 55344 21.10 99.55 62.37 14.39 93.00 7.00 1.03 0.69 65.71 4.74 28.37 1.17 25.93 17.28 25.93 79.87 1.20 17.05 0.48 97.76 58.18 2.27 11.13 4.12 4.91 1.87 6.10
66 16 63 older Male High 1 No NA BTT HMII dead Died 8 309193 61091 19.76 98.15 56.62 7.24 83.61 16.39 3.27 1.52 58.03 11.90 28.13 1.93 24.65 11.36 1.52 79.53 1.45 13.33 0.97 92.63 25.76 3.50 9.23 7.57 3.73 3.48 0.00
67 17 68 older Male Low 3 No NA DT HMII dead Died 0 147997 54649 36.93 91.75 79.63 7.62 74.86 25.14 0.60 3.45 62.45 22.87 12.30 2.38 45.98 9.28 13.61 74.33 2.67 21.19 0.84 91.31 90.15 3.58 2.25 21.79 1.36 0.58 0.34
68 17 68 older Male Low 3 No NA DT HMII dead Died 1 177687 56400 31.74 97.13 69.76 12.60 81.83 18.17 0.65 1.61 55.77 16.13 26.02 2.09 53.19 10.70 11.57 82.92 1.38 13.98 0.28 89.83 85.59 3.78 1.85 14.65 1.30 0.68 0.18
69 17 68 older Male Low 3 No NA DT HMII dead Died 3 233155 67140 28.80 98.44 74.89 5.32 81.96 18.04 1.00 1.88 49.03 16.59 32.90 1.48 49.92 7.20 12.83 83.44 1.39 12.98 0.14 89.70 69.70 4.58 1.62 13.89 1.08 1.08 0.17
70 17 68 older Male Low 3 No NA DT HMII dead Died 5 99257 29455 29.68 97.33 62.48 4.73 80.60 19.40 1.84 1.25 45.06 18.36 35.25 1.33 44.61 6.32 12.27 82.89 0.88 13.20 0.09 83.55 35.29 7.23 3.54 14.75 2.29 3.02 1.20
71 17 68 older Male Low 3 No NA DT HMII dead Died 8 76358 18702 24.49 97.32 29.48 11.10 86.29 13.71 2.43 0.50 37.92 12.13 48.32 1.63 46.45 11.70 11.35 83.96 0.20 10.50 0.06 57.82 50.00 9.95 29.60 10.20 0.64 9.36 0.80
72 18 37 younger Female High 2 Yes NA BTT PVAD dead Died s/p OHT 0 56374 12141 21.54 67.38 27.12 3.32 51.84 48.16 1.84 1.84 2.94 45.59 47.79 3.68 0.00 0.00 2.29 92.65 0.37 5.51 0.72 30.15 60.00 6.99 9.56 0.37 0.00 3.68 0.00
73 18 37 younger Female High 2 Yes NA BTT PVAD dead Died s/p OHT 3 7018 3464 49.36 89.03 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00 0.00
74 18 37 younger Female High 2 Yes NA BTT PVAD dead Died s/p OHT 8 80210 21828 27.21 85.74 20.38 3.50 61.98 38.02 8.55 2.90 5.50 34.96 56.18 3.36 0.00 0.00 3.21 89.31 0.61 2.60 1.49 23.66 0.00 21.37 9.16 1.07 0.15 5.34 1.31
75 18 37 younger Female High 2 Yes NA BTT PVAD dead Died s/p OHT 14 55068 19145 34.77 75.16 16.90 3.30 64.00 36.00 4.21 2.11 10.53 34.11 52.42 2.95 0.00 0.00 2.94 92.42 0.00 3.37 1.00 22.74 0.00 12.42 11.16 0.63 0.00 4.63 0.00
76 18 37 younger Female High 2 Yes NA BTT PVAD dead Died s/p OHT 21 76673 23615 30.80 83.39 26.27 2.06 63.55 36.45 5.42 5.42 5.17 34.48 58.13 2.22 0.00 0.00 2.72 89.41 0.00 3.45 1.95 18.23 4.55 18.23 8.87 1.23 0.25 4.68 0.00
77 19 46 younger Male High 2 No NA DT HMII alive Alive 0 226762 25868 11.41 89.41 25.31 12.05 55.06 44.94 11.31 2.69 14.79 42.57 39.95 2.69 19.87 3.88 32.88 68.95 10.48 13.57 0.59 80.62 22.67 31.91 19.38 15.58 10.37 16.55 8.18
78 19 46 younger Male High 2 No NA DT HMII alive Alive 3 244049 39386 16.14 94.94 40.46 13.23 40.93 59.07 13.36 3.17 8.33 54.90 32.28 4.49 18.52 4.21 23.83 55.18 13.87 19.39 1.44 70.53 35.67 33.62 33.47 22.86 15.26 25.81 5.11
79 19 46 younger Male High 2 No NA DT HMII alive Alive 5 344736 250121 72.55 89.83 7.54 3.81 45.84 54.16 3.96 1.53 13.39 52.70 32.18 1.72 67.14 2.92 2.44 74.68 6.81 13.91 0.21 63.33 58.02 33.05 30.94 19.83 13.96 26.74 0.22
80 19 46 younger Male High 2 No NA DT HMII alive Alive 8 216934 140093 64.58 96.37 4.78 7.22 22.49 77.51 13.49 0.64 7.22 67.61 15.17 10.00 35.88 7.82 17.22 51.28 39.40 6.56 1.33 62.59 58.06 66.67 64.59 53.00 50.26 63.78 4.97
81 19 46 younger Male High 2 No NA DT HMII alive Alive 14 215820 66785 30.94 92.54 17.20 6.01 25.88 74.12 9.15 1.45 3.09 71.72 22.44 2.74 65.43 3.42 0.04 69.38 15.93 8.85 0.74 50.58 24.07 52.78 48.88 32.28 28.19 44.04 0.04
82 19 46 younger Male High 2 No NA DT HMII alive Alive 21 670692 523014 77.98 96.05 5.07 2.64 23.61 76.39 4.11 0.49 4.79 73.53 18.68 3.00 42.40 2.72 4.49 71.95 20.81 5.51 0.45 48.22 44.62 60.13 60.26 39.07 36.82 56.38 1.49
83 20 75 older Male High 2 No NA DT HMII alive Alive 3 121794 28129 23.10 97.05 52.79 16.42 70.08 29.92 0.29 0.54 50.71 21.46 22.96 4.86 34.18 12.51 23.49 39.78 0.31 55.18 0.35 96.01 16.67 1.74 5.27 23.38 2.95 1.00 0.10
84 20 75 older Male High 2 No NA DT HMII alive Alive 5 82157 17295 21.05 99.20 72.54 7.41 68.21 31.79 1.10 0.55 48.47 22.90 23.60 5.04 31.90 13.57 18.10 37.84 0.47 56.73 0.35 92.37 14.29 2.05 7.79 26.20 4.72 1.42 1.02
85 20 75 older Male High 2 No NA DT HMII alive Alive 8 154810 27067 17.48 98.04 72.18 5.03 74.91 25.09 0.90 1.65 58.05 14.76 20.82 6.37 23.43 25.71 22.57 44.94 0.60 49.36 0.71 86.59 22.73 3.30 7.64 20.15 4.87 1.72 0.00
86 20 75 older Male High 2 No NA DT HMII alive Alive 21 197936 15446 7.80 95.86 50.35 5.44 79.13 20.87 3.98 3.11 52.30 15.90 29.81 1.99 20.23 5.78 24.86 40.62 1.12 50.06 1.89 62.36 8.00 7.20 9.94 13.17 3.60 5.09 1.55
87 21 30 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 0 114586 35996 31.41 100.00 59.03 11.31 75.87 24.13 1.52 7.69 69.83 11.67 10.20 8.30 27.80 29.94 2.04 84.64 0.64 10.54 0.61 91.57 61.34 4.99 7.15 11.57 5.92 2.53 1.92
88 21 30 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 3 349054 129419 37.08 95.25 62.35 8.81 79.64 20.36 2.72 4.45 57.70 13.93 25.89 2.48 57.14 13.49 8.69 74.29 3.76 12.35 1.41 85.87 62.53 3.41 5.43 7.67 1.40 1.92 1.06
89 21 30 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 5 447309 176789 39.52 93.00 53.57 9.32 80.92 19.08 2.16 5.73 56.72 12.51 28.43 2.34 47.68 12.19 11.00 73.80 4.44 12.53 2.57 83.88 49.09 5.24 6.82 8.32 2.18 2.28 1.72
90 21 30 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 8 558290 301468 54.00 95.62 47.73 7.04 79.13 20.87 3.43 7.49 53.19 13.48 30.32 3.01 36.36 14.54 12.10 70.55 5.15 11.89 3.51 79.79 52.89 6.98 8.50 7.82 1.92 3.09 1.96
91 21 30 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 14 479402 362845 75.69 97.95 44.30 6.30 76.53 23.47 3.90 11.10 34.99 12.04 50.31 2.66 21.93 14.88 10.05 70.21 4.66 7.31 7.73 49.96 29.38 13.44 17.37 5.64 3.01 6.95 2.28
92 21 30 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 21 334942 191418 57.15 93.28 70.35 4.04 70.89 29.11 6.20 17.70 56.45 19.57 20.44 3.53 38.91 14.64 10.35 55.11 6.61 14.34 8.70 80.96 58.92 7.58 7.96 10.11 2.08 2.08 2.48
93 22 67 older Male High 2 Yes NA BTT HMII dead Died 0 129980 41600 32.00 91.40 39.03 2.85 49.63 50.37 7.02 11.55 28.37 29.11 20.79 21.72 19.43 28.34 9.80 69.59 4.90 20.33 6.13 31.15 33.60 20.15 22.27 0.28 4.25 16.54 2.09
94 22 67 older Male High 2 Yes NA BTT HMII dead Died 1 59478 20182 33.93 89.07 44.50 10.08 60.93 39.07 2.10 2.70 42.88 13.69 17.66 25.77 13.06 46.12 3.13 81.35 1.16 16.11 0.83 37.58 6.12 5.52 11.64 0.06 1.49 2.76 4.07
95 22 67 older Male High 2 Yes NA BTT HMII dead Died 5 26646 9754 36.61 98.69 58.32 7.57 64.88 35.12 2.06 0.55 46.50 12.62 17.70 23.18 15.07 53.31 12.50 83.40 1.10 15.64 1.76 28.12 0.00 7.54 16.46 0.55 1.92 6.04 13.00
96 22 67 older Male High 2 Yes NA BTT HMII dead Died 8 73641 21944 29.80 96.79 48.91 4.47 62.74 37.26 2.84 3.16 37.47 17.79 24.95 19.79 19.46 38.11 9.46 79.79 1.89 13.58 1.24 27.05 20.00 14.11 13.37 0.21 2.74 10.53 3.31
97 23 68 older Male High 2 No No BTT HMII alive Alive s/p OHT 0 44004 20612 46.84 43.65 15.61 21.54 49.64 50.36 10.37 0.57 15.94 45.15 33.18 5.73 24.69 5.68 20.25 44.58 1.55 36.17 0.88 67.23 27.27 32.35 15.48 20.02 4.49 10.94 2.76
98 23 68 older Male High 2 No No BTT HMII alive Alive s/p OHT 1 99731 23424 23.49 91.81 41.33 10.89 35.24 64.76 11.83 1.62 18.71 59.21 16.19 5.89 28.16 5.89 18.53 39.98 4.66 42.16 2.33 84.07 2.63 22.55 15.21 40.28 6.07 9.70 3.96
99 23 68 older Male High 2 No No BTT HMII alive Alive s/p OHT 8 129637 50001 38.57 90.65 44.41 6.05 30.66 69.34 17.12 4.05 10.62 66.50 19.60 3.28 26.41 3.46 13.37 39.74 7.15 39.49 11.91 71.97 2.70 32.48 33.14 43.91 7.81 22.70 3.82
100 23 68 older Male High 2 No No BTT HMII alive Alive s/p OHT 14 134368 46363 34.50 88.32 34.25 4.63 32.47 67.53 12.86 8.38 8.91 64.89 22.88 3.32 24.78 3.89 13.58 36.16 6.01 40.33 13.77 72.11 1.89 29.89 26.36 42.80 6.17 15.97 3.29
101 24 64 older Male High 2 No Yes BTT HMII alive Alive s/p OHT 0 59881 8971 14.98 97.40 50.06 9.57 75.84 24.16 3.59 8.13 64.83 13.40 15.67 6.10 19.31 18.81 16.83 76.67 6.10 6.82 2.00 88.16 47.06 10.53 11.24 10.41 9.93 6.94 11.93
102 24 64 older Male High 2 No Yes BTT HMII alive Alive s/p OHT 1 60193 23874 39.66 96.72 58.83 16.24 80.37 19.63 1.87 3.31 68.58 10.62 14.72 6.08 21.33 22.83 3.80 80.61 3.63 8.03 0.33 94.16 57.26 2.96 5.01 10.19 4.72 1.41 1.24
103 24 64 older Male High 2 No Yes BTT HMII alive Alive s/p OHT 5 135335 27834 20.57 97.35 44.75 5.61 63.22 36.78 8.49 17.63 39.01 23.16 30.72 7.11 17.53 10.38 12.16 42.70 4.93 19.28 4.68 67.37 4.10 26.64 17.76 19.21 14.54 14.47 10.60
104 24 64 older Male High 2 No Yes BTT HMII alive Alive s/p OHT 8 247113 62042 25.11 96.41 43.93 14.82 78.42 21.58 2.52 5.24 65.35 12.22 17.53 4.90 20.39 17.20 2.93 77.86 2.11 8.78 2.94 84.64 9.70 10.64 13.93 9.90 11.59 8.48 3.51
105 24 64 older Male High 2 No Yes BTT HMII alive Alive s/p OHT 14 142053 56983 40.11 93.01 62.72 12.87 85.10 14.90 0.69 1.64 68.02 8.40 21.03 2.55 17.03 13.68 1.08 83.33 0.37 7.38 0.26 79.00 16.07 3.01 5.85 5.16 5.24 1.77 0.70
106 24 64 older Male High 2 No Yes BTT HMII alive Alive s/p OHT 21 60384 17029 28.20 93.65 40.96 21.36 83.80 16.20 1.47 2.32 63.52 9.69 23.80 2.99 23.55 15.58 0.91 78.16 1.09 8.13 1.00 83.65 17.72 4.70 8.42 5.99 7.22 3.76 1.56
107 25 25 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 0 47645 37377 78.45 94.61 2.47 9.41 78.56 17.02 1.71 4.63 34.94 16.09 44.77 4.21 4.35 28.99 59.42 67.26 3.10 26.22 7.50 72.82 53.90 13.50 28.20 3.04 7.70 7.67 0.96
108 25 25 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 1 87452 43025 49.20 91.54 1.98 26.71 92.06 6.76 0.77 4.67 47.21 4.70 45.17 2.92 6.06 18.18 60.61 81.01 3.95 14.10 12.04 85.97 79.84 5.82 14.15 1.11 4.77 2.84 0.00
109 25 25 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 3 124056 80873 65.19 95.41 0.37 14.36 90.47 5.98 0.39 0.60 41.84 4.26 51.08 2.81 18.82 35.29 43.53 86.97 0.10 11.32 1.10 75.49 12.12 3.80 14.94 0.24 0.88 3.25 0.00
110 25 25 younger Female Low 3 No Yes BTT HMII alive Alive s/p OHT 8 87840 57990 66.02 95.60 0.11 15.80 89.00 9.49 1.56 5.02 45.89 7.18 43.54 3.39 6.10 19.51 56.71 84.23 3.14 9.95 10.70 72.13 52.73 7.06 17.53 0.74 3.13 5.13 0.13
111 26 49 younger Male Low 3 No NA DT HMII alive Alive 0 196315 122053 62.17 95.37 60.66 22.66 91.10 8.90 0.45 4.95 85.57 6.30 5.42 2.71 13.13 28.32 20.05 86.32 2.45 5.42 0.96 98.68 94.19 2.50 2.96 4.63 1.24 0.33 0.89
112 26 49 younger Male Low 3 No NA DT HMII alive Alive 1 91909 37622 40.93 96.10 43.22 38.85 95.48 4.52 0.14 2.19 85.07 3.09 10.33 1.51 21.03 27.23 9.08 88.84 1.13 3.91 0.23 98.77 88.27 1.16 1.16 2.34 0.38 0.21 0.45
113 26 49 younger Male Low 3 No NA DT HMII alive Alive 5 153762 86281 56.11 95.83 65.54 24.24 91.12 8.88 0.65 2.84 81.43 6.30 9.62 2.65 21.14 28.73 15.89 83.00 1.97 6.74 0.36 98.22 93.16 3.20 2.29 5.47 1.11 0.77 1.56
114 26 49 younger Male Low 3 No NA DT HMII alive Alive 8 82406 47485 57.62 92.71 65.75 17.46 91.79 8.21 1.08 1.99 76.35 5.70 15.41 2.54 20.56 26.95 11.99 78.67 1.35 7.97 0.28 97.55 85.62 2.65 2.56 5.19 1.18 0.86 1.79
115 27 62 older Male High 2 No NA BTT CMAG dead Died 0 105070 38521 36.66 98.39 47.84 12.88 60.82 39.18 11.08 2.89 33.87 29.43 26.68 10.01 19.72 21.92 2.99 59.35 14.64 20.23 3.45 81.40 46.10 24.47 24.27 24.04 11.00 17.02 1.51
116 27 62 older Male High 2 No NA BTT CMAG dead Died 1 162336 73080 45.02 99.08 49.35 13.95 74.24 25.76 7.77 3.84 45.21 19.33 28.85 6.61 19.16 23.45 0.81 65.34 14.40 13.50 3.17 86.68 71.39 23.26 17.70 15.21 6.93 13.29 0.88
117 27 62 older Male High 2 No NA BTT CMAG dead Died 5 111264 60487 54.36 99.15 45.94 16.13 73.16 26.84 4.64 1.82 52.04 20.07 20.94 6.95 22.93 23.47 0.50 77.11 9.30 9.58 1.20 88.92 81.82 22.20 19.01 16.98 6.89 13.57 1.36
118 27 62 older Male High 2 No NA BTT CMAG dead Died 8 165870 67935 40.96 98.88 43.01 15.86 70.23 29.77 7.85 1.88 40.39 20.49 29.66 9.46 22.70 29.00 1.88 68.61 9.60 16.11 1.47 89.57 69.00 20.12 16.59 17.10 7.02 12.40 2.06
119 28 72 older Male Low 4 Yes NA BTT HMII dead Died 3 424375 21971 5.18 71.78 56.98 18.53 79.84 20.16 3.08 0.51 57.91 14.07 24.06 3.97 31.74 15.95 13.65 57.94 0.65 28.51 0.17 93.57 53.33 2.87 9.69 15.61 3.87 2.16 3.99
120 28 72 older Male Low 4 Yes NA BTT HMII dead Died 8 177249 103442 58.36 76.40 15.43 5.64 65.13 34.87 5.21 0.20 35.56 27.91 31.90 4.62 26.60 8.28 13.24 67.89 0.29 20.13 0.00 78.60 22.22 5.50 16.09 25.56 6.80 4.73 0.89

Identification of clusters

We double standardized the data and biclustered it. We found 3 clusters of patients, and 3 clusters of B-cells.

require(pheatmap, quietly = T)
require(RColorBrewer, quietly = T)


colorfun <- function(nlevels, ...){
    if(nlevels > 10) return(standardColors(nlevels))
    if(nlevels <=10) return(pal_d3(...)(nlevels))
} 

annotation.col <- df.raw[,c(1,3,4,6,7,8,9,10,11,12)]
annotation.colors <- lapply(colnames(annotation.col), function(nam){
    colrs <- colorfun(nlevels(annotation.col[[nam]]))
    names(colrs) <- levels(annotation.col[[nam]])
    return(colrs)
})
names(annotation.colors) <- colnames(annotation.col)

df.patient <- double_standardize(t(df.raw[,c(2,14:42)]))
pheatmap(df.patient, 
         color = colorRampPalette(rev(brewer.pal(n = 7, name = "RdBu")))(100),
         breaks = c(seq(-max(abs(df.patient), na.rm=T), 0, length.out = 50), seq(0.001, max(abs(df.patient), na.rm=T), length.out = 50)),
         #scale = "row",
         fontsize = 7,
         cutree_rows = 3,
         cutree_cols = 3,
         annotation_col = annotation.col,
         annotation_colors = annotation.colors
)

Metadata associations

We computed the statistical associations between sample groups specified by the metadata.

Survival

We computed metadata factors that were statistically associated with survival.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("Survival","Outcome", "LowIntermacs"))], 
                      function(this_factor){
                          table(df[,this_factor], df$Survival)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 5.200 2.650 10.40 1.93e-07 1.11e-06 RVAD
2 NA NA NA 3.16e-07 1.11e-06 Device Type
3 NA NA NA 0.00153 0.00358 InterMACS
4 2.770 0.931 9.01 0.0513 0.0899 Sensitized
5 1.500 0.797 2.84 0.232 0.273 AgeGreater60
6 0.668 0.328 1.37 0.234 0.273 Sex
7 0.876 0.392 1.89 0.856 0.856 VAD Indication

InterMACS

We computed metadata factors that were statistically associated with a binarized InterMACS score.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("LowIntermacs", "InterMACS"))], 
                      function(this_factor){
                          table(df[,this_factor], df$LowIntermacs)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 NA NA NA 2.4e-13 1.92e-12 Device Type
2 NA NA NA 2.99e-05 0.000119 Outcome
3 3.190 1.690 6.10 0.00013 0.000346 AgeGreater60
4 1.860 0.965 3.67 0.0493 0.0986 RVAD
5 0.573 0.271 1.21 0.113 0.181 VAD Indication
6 1.330 0.575 3.10 0.557 0.701 Sensitized
7 1.220 0.596 2.47 0.613 0.701 Sex
8 0.955 0.504 1.82 0.88 0.88 Survival

Sensitization

We computed metadata factors that were statistically associated with sensitization.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("Sensitized"))], 
                      function(this_factor){
                          table(df[,this_factor], df$Sensitized)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 NA NA NA 1.1e-06 5.83e-06 Outcome
2 0.110 0.0347 0.307 1.3e-06 5.83e-06 Sex
3 NA NA NA 0.0129 0.0388 InterMACS
4 NA NA NA 0.019 0.0427 Device Type
5 0.404 0.1630 0.967 0.0295 0.0531 RVAD
6 2.770 0.9310 9.010 0.0513 0.077 Survival
7 1.330 0.5750 3.100 0.557 0.716 LowIntermacs
8 1.250 0.5310 2.940 0.69 0.776 AgeGreater60
9 0.000 0.0000 Inf 1 1 VAD Indication

Sex

We computed metadata factors that were statistically associated with sex.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("Sex"))], 
                      function(this_factor){
                          table(df[,this_factor], df$Sex)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 0.110 0.0347 0.307 1.3e-06 1.17e-05 Sensitized
2 NA NA NA 2.89e-06 1.3e-05 Device Type
3 NA NA NA 2.21e-05 6.64e-05 Outcome
4 4.030 1.9200 8.870 6.13e-05 0.000138 AgeGreater60
5 NA NA NA 0.0309 0.0555 InterMACS
6 1.870 0.7420 5.380 0.227 0.263 VAD Indication
7 0.668 0.3280 1.370 0.234 0.263 RVAD
8 0.668 0.3280 1.370 0.234 0.263 Survival
9 1.220 0.5960 2.470 0.613 0.613 LowIntermacs

Age

We computed metadata factors that were statistically associated with age.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("AgeGreater60"))], 
                      function(this_factor){
                          table(df[,this_factor], df$AgeGreater60)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 NA NA NA 9.53e-11 8.58e-10 Outcome
2 NA NA NA 5.7e-09 2.56e-08 InterMACS
3 NA NA NA 2.35e-08 7.04e-08 Device Type
4 0.216 0.108 0.421 1.25e-06 2.82e-06 RVAD
5 4.030 1.920 8.870 6.13e-05 0.00011 Sex
6 3.190 1.690 6.100 0.00013 0.000195 LowIntermacs
7 1.990 0.931 4.430 0.0577 0.0742 VAD Indication
8 1.500 0.797 2.840 0.232 0.261 Survival
9 1.250 0.531 2.940 0.69 0.69 Sensitized

VAD Indication

We computed metadata factors that were statistically associated with VAD Indication.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("VAD Indication"))], 
                      function(this_factor){
                          table(df[,this_factor], df$`VAD Indication`)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = F, B = 10000)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 NA NA NA 4.18e-30 3.76e-29 Outcome
2 0.000 0.000 0.116 5.89e-10 2.65e-09 RVAD
3 NA NA NA 6.72e-05 0.000201 Device Type
4 NA NA NA 0.0272 0.0612 InterMACS
5 1.990 0.931 4.430 0.0577 0.104 AgeGreater60
6 0.573 0.271 1.210 0.113 0.17 LowIntermacs
7 1.870 0.742 5.380 0.227 0.292 Sex
8 0.876 0.392 1.890 0.856 0.963 Survival
9 0.000 0.000 Inf 1 1 Sensitized

RVAD

We computed metadata factors that were statistically associated with RVAD.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("RVAD"))], 
                      function(this_factor){
                          table(df[,this_factor], df$RVAD)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 NA NA NA 6.86e-20 6.17e-19 Device Type
2 NA NA NA 1.39e-13 6.27e-13 Outcome
3 0.000 0.000 0.116 5.89e-10 1.77e-09 VAD Indication
4 5.200 2.650 10.400 1.93e-07 4.35e-07 Survival
5 0.216 0.108 0.421 1.25e-06 2.25e-06 AgeGreater60
6 NA NA NA 2.8e-05 4.21e-05 InterMACS
7 0.404 0.163 0.967 0.0295 0.0379 Sensitized
8 1.860 0.965 3.670 0.0493 0.0555 LowIntermacs
9 0.668 0.328 1.370 0.234 0.234 Sex

Device-type

We computed metadata factors that were statistically associated with Device Type.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[!(names(isfactor) %in% c("Device Type"))], 
                      function(this_factor){
                          table(df[,this_factor], df$`Device Type`)
                      })

fisher.p <- sapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = T, B = 10000)
    ft$p
})

fisher.OR.0 <- lapply(contingency, function(this_table){
    ft <- fisher.test(this_table, simulate.p.value = T, B = 10000)
    if(is.null(ft$estimate)) {
        ft$estimate <- NA
        ft$conf.int <- c(NA,NA)
    }
    ret <- list(OR = ft$estimate, lowerCL = ft$conf.int[1], upperCL = ft$conf.int[2])
    return(ret)
})
fisher.OR <- as.data.frame(do.call(rbind, fisher.OR.0))
fisher.OR$pvalue <- fisher.p
fisher.OR$qvalue <- p.adjust(fisher.p, method = "BH")
fisher.all <- as.data.frame(sapply(fisher.OR, as.numeric))
rownames(fisher.all) <- rownames(fisher.OR)

qtable <- signif(fisher.all[order(fisher.all$pvalue),],3)

qtable %>%
    mutate(
        `Factor` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>% 
    kable(escape = F, row.names = T) %>%
    kable_styling(bootstrap_options = c("striped", 
                                        "hover", 
                                        "condensed",
                                        "responsive"),
                  font_size = 12)
OR lowerCL upperCL pvalue qvalue Factor
1 NA NA NA 1e-04 0.000129 AgeGreater60
2 NA NA NA 1e-04 0.000129 Sex
3 NA NA NA 1e-04 0.000129 LowIntermacs
4 NA NA NA 1e-04 0.000129 InterMACS
5 NA NA NA 1e-04 0.000129 RVAD
6 NA NA NA 1e-04 0.000129 Survival
7 NA NA NA 1e-04 0.000129 Outcome
8 NA NA NA 2e-04 0.000225 VAD Indication
9 NA NA NA 0.0174 0.0174 Sensitized

Visualization of metadata dependencies

We used the GGally package to visualize pairwise dependencies between sample groups, separating out survivors and non-survivors (indicated by color).

df.temp <- df
colnames(df.temp) <- make.names(colnames(df), unique = T)
#invisible(
suppressWarnings(
    suppressMessages(
        ggpairs(df.temp, 
                mapping = aes(color = Survival), 
                columns = colnames(df.temp)[c(2,4,6,7,8,9,10,11)]) 
    ))

#)

Variability of MCS devices

PCA

Using Principal Component Analysis (PCA), we saw large variability between device types, compared to the variability within device types. We also saw large variability between individual patients. None of the other features were clearly separable.

suppressMessages(require(ggbiplot, quietly = T))
require(ggsci, quietly = T)

isna <- unique(unlist(apply(df[,14:42], 2, function(x) which(is.na(x)))))
pca <- prcomp(double_standardize(df[-isna, 14:42]), center = TRUE, scale. = TRUE)

colorfun <- function(grouping, ...){
    if(nlevels(factor(grouping)) > 10) scale_color_discrete(...)
    if(nlevels(factor(grouping)) <=10) scale_color_d3(...)
} 

plotvars <- names(df)[c(10,1:9,11:17)]
plots.pca <- mclapply(plotvars, function(this_var){
    this_groups <- df[-isna, this_var]
    if(is.numeric(this_groups)) this_groups <- cut(this_groups, 4)
    
    ggbiplot(pca,
             groups = this_groups,
             ellipse = TRUE,
             alpha = 0.3,
             varname.size = 1.2) + 
        colorfun(this_groups, name = this_var) +
        ggtitle(paste0(this_var, " variability")) +
        theme_classic() 
}, mc.cores = detectCores()-1)
names(plots.pca) <- plotvars
#plots.pca$PatientID
#plots.pca$`Device Type`

for(ii in 1:length(plots.pca)){
    cat("  \n###", names(plots.pca)[ii], "\n")
    suppressWarnings(print(plots.pca[[ii]]))
    cat("  \n")
}

Device Type

PatientID

Age

AgeGreater60

Sex

LowIntermacs

InterMACS

RVAD

Sensitized

VAD Indication

Survival

Outcome

Time

num Total PBMC

num lymph

lymph

live lymph

Fisher’s exact test

We analyzed the dependency of device types on the other discrete variables using Fisher’s exact test. We found several dependencies (Benjamini-Hochberg \(q<0.05\)), which we have shown below.

isfactor <- which(sapply(df, is.factor))[-1]
contingency <- lapply(isfactor[names(isfactor) != "Device Type"], 
                      function(this_factor){
                          table(df[,this_factor], df$`Device Type`)
                      })

fisher.p <- sapply(contingency, function(this_table){
    fisher.test(this_table, simulate.p.value = T, B = 10000)$p
})
fisher.q <- p.adjust(fisher.p, method = "BH")

signif.ix <- which(fisher.q < 0.05)
signif.order <- sort(fisher.q[signif.ix], index.return = T)$ix
for(this_ix in signif.ix[signif.order]){
    cat("  \n###", names(contingency)[this_ix], "\n")
    cat(paste0("Benjamini-Hochberg qvalue = ", signif(fisher.q[this_ix], 2)),
        ". \n")
    print(kable(contingency[this_ix][[1]],  row.names = T) %>%
              kable_styling(bootstrap_options = c("striped", 
                                                  "hover", 
                                                  "condensed",
                                                  "responsive"),
                            font_size = 12)
    )
    cat("  \n")
}

AgeGreater60

Benjamini-Hochberg qvalue = 0.00013 .
HMII CMAG HVAD PVAD TAH
younger 56 0 7 21 7
older 84 7 7 0 7

Sex

Benjamini-Hochberg qvalue = 0.00013 .
HMII CMAG HVAD PVAD TAH
Female 35 0 0 14 0
Male 105 7 14 7 14

LowIntermacs

Benjamini-Hochberg qvalue = 0.00013 .
HMII CMAG HVAD PVAD TAH
Low 63 0 14 0 0
High 77 7 0 21 14

InterMACS

Benjamini-Hochberg qvalue = 0.00013 .
HMII CMAG HVAD PVAD TAH
1 7 0 0 0 7
2 70 7 0 21 7
3 56 0 14 0 0
4 7 0 0 0 0

RVAD

Benjamini-Hochberg qvalue = 0.00013 .
HMII CMAG HVAD PVAD TAH
No 105 7 14 0 0
Yes 35 0 0 21 14

Survival

Benjamini-Hochberg qvalue = 0.00013 .
HMII CMAG HVAD PVAD TAH
alive 98 0 14 7 7
dead 42 7 0 14 7

Outcome

Benjamini-Hochberg qvalue = 0.00013 .
HMII CMAG HVAD PVAD TAH
Alive 28 0 0 0 0
Alive s/p OHT 70 0 14 7 7
Died 35 7 0 7 0
Died post OHT 7 0 0 0 0
Died s/p OHT 0 0 0 7 7

VAD Indication

Benjamini-Hochberg qvalue = 0.00045 .
HMII CMAG HVAD PVAD TAH
BTT 98 7 14 21 14
DT 42 0 0 0 0

Sensitized

Benjamini-Hochberg qvalue = 0.018 .
HMII CMAG HVAD PVAD TAH
No 42 0 7 0 7
Yes 28 0 7 7 7

One-way repeated measures anova

We attempted to analyze the differences in B-cell levels across device types using a linear mixed effect model. Here we report variables that had a statistically significant variance (Benjamini-Hochberg \(q<0.05\)) across devices, time, or their interaction.

suppressMessages(require(lmerTest, quietly = T))
suppressMessages(require(car, quietly = TRUE))

df.device <- df
colnames(df.device) <- make.names(colnames(df), unique = T)

varnames <- colnames(df.device)[14:42]
models.device <- mclapply(varnames, function(this_var){
    this_formula <- as.formula(paste0(this_var, " ~ Device.Type + (1|PatientID)"))
    invisible(suppressMessages(this_model <- lmer(this_formula, data = df.device)))
    this_anova <- Anova(this_model, type = 2)
    pvals <- this_anova$`Pr(>Chisq)`[c(1)]
    return(list(model = this_model,
                pvals = pvals))
}, mc.cores = detectCores()-1)
names(models.device) <- colnames(df)[14:42]

pvals <- do.call(rbind, lapply(models.device, function(x) x$pvals))
rownames(pvals) <- colnames(df)[14:42]
colnames(pvals) <- c("pvalue")

qBH <- matrix(p.adjust(pvals, method = "BH"), nrow = nrow(pvals))
rownames(qBH) <- colnames(df)[14:42]
colnames(qBH) <- c("qvalue")

sigvars <- apply(apply(pvals, 2, function(x) x<=0.05), 1, function(x) sum(x) > 0)

sig.qtable <- cbind(pvals,qBH)[sigvars,,drop=F][order(apply(pvals[sigvars,,drop=F], 1, min)),,drop=F]

qtable <- as.data.frame(signif(sig.qtable, 2))
qtable %>%
    mutate(
        `B-cell` = row.names(.),
        pvalue = cell_spec(pvalue, color = ifelse(qtable$pvalue > 0.05, "grey", "red")),
        qvalue = cell_spec(qvalue, color = ifelse(qtable$qvalue > 0.05, "grey", "red"))
    ) %>%
    kable( escape = F,
           digits = 20,
           row.names = T,
           caption = "Significant device-type q-values") %>%
    kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
                  font_size = 12) %>%
    scroll_box(width = "100%")
Significant device-type q-values
pvalue qvalue B-cell
1 0.011 0.31 CD27-IgD+ mature naive
2 0.045 0.45 CD19+CD5+

We plotted the average across time for each of the B-cells that showed a statistically significant effect across devices in the above mixed effect models. We drew attention to specific features that induced the positive test result, by listing the model parameters with \(p<0.05\) in the fit. Note that the reference level for the devices is the HeartMate-II (HMII).

require(reshape2, quietly = T)
df.long <- melt(df, id.vars = colnames(df)[1:13])
groups <- make.names(c("Device Type"))
names(df.long) <- make.names(names(df.long))

invisible(suppressMessages(require(Hmisc, quietly = T)))
stat_sum_df <- function(fun, geom="errorbar", ...) {
    stat_summary(fun.data = fun, geom = geom, width = 1, ...)
}

plots.ts <- mclapply(rownames(qtable), function(this_var){
    lapply(groups, function(this_groups){
        ggplot(subset(df.long, df.long$variable == this_var)) +
            aes(x = Time, y = value, group = PatientID) +
            aes_string(color = this_groups, fill = this_groups) +
            geom_line(alpha = 0) + 
            geom_point(alpha = 0) + 
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("point"), position = position_dodge(.5)) + 
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("line"), position = position_dodge(.5)) + 
            stat_sum_df(function(x) mean_cl_normal(x, conf.int = 0.68), mapping = aes_string(group = this_groups), position = position_dodge(.5)) + 
            #stat_smooth(aes_string(group = this_groups), method = "loess", span = 1) +
            scale_color_d3() + scale_fill_d3() +
            xlab("Time (days after surgery)") +
            ylab(this_var) + 
            ggtitle(paste(this_var)) +
            theme_classic()
    })
}, mc.cores = detectCores()-1)
names(plots.ts) <- rownames(qtable)

for(ii in c(1:length(plots.ts))){
    for(jj in 1:length(plots.ts[[ii]])){
        sumtable <- suppressMessages(summary(models.device[[rownames(qtable)[ii]]]$model))
        sumtable <- sumtable$coefficients[-1, c(1:5)] # drop intercept
        sigsum <- sumtable[sumtable[,5] <= 0.05, , drop = F]
        cat("  \n###", rownames(qtable)[ii], "\n")
        print(kable(sigsum[order(sigsum[,5]),,drop=F], row.names = T) %>%
                  kable_styling(bootstrap_options = c("striped", 
                                                      "hover", 
                                                      "condensed",
                                                      "responsive"),
                                font_size = 12)
        )
        cat("  \n")
        suppressWarnings(print(plots.ts[[ii]][[jj]]))
        cat("  \n")
    }
}

CD27-IgD+ mature naive

Estimate Std. Error df t value Pr(>|t|)
Device.TypePVAD -33.85486 10.06687 23.54824 -3.362998 0.0026287

CD19+CD5+

Estimate Std. Error df t value Pr(>|t|)
Device.TypePVAD 23.16366 8.615765 23.18582 2.688521 0.0130627

Mixed-effect model

We attempted to analyze the differences in B-cell levels across device types using a linear mixed effect model. We used time as a continuous variable, and included the interaction term. Here we report variables that had a statistically significant variance (Benjamini-Hochberg \(q<0.05\)) across devices, time, or their interaction.

suppressMessages(require(lmerTest, quietly = T))
suppressMessages(require(car, quietly = TRUE))

df.device <- df
colnames(df.device) <- make.names(colnames(df), unique = T)

varnames <- colnames(df.device)[14:42]
models.device <- mclapply(varnames, function(this_var){
    this_formula <- as.formula(paste0(this_var, " ~ Device.Type * Time + (1|PatientID)"))
    invisible(suppressMessages(this_model <- lmer(this_formula, data = df.device)))
    this_anova <- Anova(this_model, type = 2)
    pvals <- this_anova$`Pr(>Chisq)`[c(1,2,3)]
    return(list(model = this_model,
                pvals = pvals))
}, mc.cores = detectCores()-1)
names(models.device) <- colnames(df)[14:42]

pvals <- do.call(rbind, lapply(models.device, function(x) x$pvals))
rownames(pvals) <- varnames
colnames(pvals) <- c("Device", "Time", "Device:Time")

qBH <- matrix(p.adjust(pvals, method = "BH"), nrow = nrow(pvals))
rownames(qBH) <- colnames(df)[14:42]
colnames(qBH) <- colnames(pvals)

sigvars <- apply(apply(qBH, 2, function(x) x<=0.05), 1, function(x) sum(x) > 0)
sig.qtable <- qBH[sigvars,][order(apply(qBH[sigvars,], 1, min)),]

qtable <- as.data.frame(signif(sig.qtable, 2))
qtable %>%
    mutate(
        `B-cell` = row.names(.),
        Device = cell_spec(Device, color = ifelse(qtable$Device > 0.05, "grey", "red")),
        `Device:Time` = cell_spec(`Device:Time`, color = ifelse(qtable$`Device:Time` > 0.05, "grey", "red")),
        `Time` = cell_spec(`Time`, color = ifelse(qtable$`Time` > 0.05, "grey", "red"))
    ) %>%
    kable( escape = F,
           digits = 20,
           row.names = T,
           caption = "Significant device-type q-values") %>%
    kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
                  font_size = 12) %>%
    scroll_box(width = "100%")
Significant device-type q-values
Device Time Device:Time B-cell
1 0.6 0.0023 0.007 CD19+CD11b+
2 0.5 0.025 0.82 CD268 of +27-38++transitional
3 0.82 0.025 0.82 CD19+CD5+CD11b+

We plotted the average across time for each of the B-cells that showed a statistically significant effect across devices in the above mixed effect models. We drew attention to specific features that induced the positive test result, by listing the model parameters with \(p<0.05\) in the multivariate fit. Note that the reference level for the devices is the HeartMate-II (HMII).

require(reshape2, quietly = T)
df.long <- melt(df, id.vars = colnames(df)[1:13])
groups <- make.names(c("Device Type"))
names(df.long) <- make.names(names(df.long))

invisible(suppressMessages(require(Hmisc, quietly = T)))
stat_sum_df <- function(fun, geom="errorbar", ...) {
    stat_summary(fun.data = fun, geom = geom, width = 1, ...)
}

plots.ts <- mclapply(rownames(qtable), function(this_var){
    lapply(groups, function(this_groups){
        ggplot(subset(df.long, df.long$variable == this_var)) +
            aes(x = Time, y = value, group = PatientID) +
            aes_string(color = this_groups, fill = this_groups) +
            geom_line(alpha = 0) + 
            geom_point(alpha = 0) + 
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("point"), position = position_dodge(.5)) + 
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("line"), position = position_dodge(.5)) + 
            stat_sum_df(function(x) mean_cl_normal(x, conf.int = 0.68), mapping = aes_string(group = this_groups), position = position_dodge(.5)) + 
            #stat_smooth(aes_string(group = this_groups), method = "loess", span = 1) +
            scale_color_d3() + scale_fill_d3() +
            xlab("Time (days after surgery)") +
            ylab(this_var) + 
            ggtitle(paste(this_var)) +
            theme_classic()
    })
}, mc.cores = detectCores()-1)
names(plots.ts) <- rownames(qtable)

for(ii in c(1:length(plots.ts))){
    for(jj in 1:length(plots.ts[[ii]])){
        sumtable <- suppressMessages(summary(models.device[[rownames(qtable)[ii]]]$model))
        sumtable <- sumtable$coefficients[-1, c(1:5)] # drop intercept
        sigsum <- sumtable[sumtable[,5] <= 0.05, , drop = F]
        cat("  \n###", rownames(qtable)[ii], "\n")
        print(kable(sigsum[order(sigsum[,5]),,drop=F], row.names = T) %>%
                  kable_styling(bootstrap_options = c("striped", 
                                                      "hover", 
                                                      "condensed",
                                                      "responsive"),
                                font_size = 12)
        )
        cat("  \n")
        suppressWarnings(print(plots.ts[[ii]][[jj]]))
        cat("  \n")
    }
}

CD19+CD11b+

Estimate Std. Error df t value Pr(>|t|)
Device.TypeTAH:Time 2.1101624 0.5538816 89.91928 3.809771 0.0002541
Device.TypePVAD:Time 0.8054722 0.2943289 88.29844 2.736640 0.0075041
Time 0.2628129 0.1188074 89.80591 2.212092 0.0294973

CD268 of +27-38++transitional

Estimate Std. Error df t value Pr(>|t|)
Time -0.7873216 0.3030382 89.7994 -2.598094 0.0109556

CD19+CD5+CD11b+

Estimate Std. Error df t value Pr(>|t|)
Time 0.290629 0.1150892 91.27507 2.525249 0.0132869

Two-way repeated measures anova

We attempted to analyze the differences in B-cell levels across device types using a mixed effect model (a.k.a. two-way repeated measures ANOVA). Here we report variables that had a statistically significant variance (Benjamini-Hochberg \(q<0.05\)) across devices, times, or their interaction. As there were 5 device types and 7 timepoints, but only 120 samples, this model is severely underpowered. The posterior belief in any of these results should therefore be quite small (as a consequence of Bayes rule).

suppressMessages(require(lmerTest, quietly = T))
suppressMessages(require(car, quietly = TRUE))

df.device <- df
colnames(df.device) <- make.names(colnames(df), unique = T)

varnames <- colnames(df.device)[14:42]
models.device <- mclapply(varnames, function(this_var){
    this_formula <- as.formula(paste0(this_var, " ~ Device.Type * factor(Time) + (1|PatientID)"))
    invisible(suppressMessages(this_model <- lmer(this_formula, data = df.device)))
    this_anova <- Anova(this_model, type = 2)
    pvals <- this_anova$`Pr(>Chisq)`[c(1,2,3)]
    return(list(model = this_model,
                pvals = pvals))
}, mc.cores = detectCores()-1)
names(models.device) <- colnames(df)[14:42]

pvals <- do.call(rbind, lapply(models.device, function(x) x$pvals))
rownames(pvals) <- varnames
colnames(pvals) <- c("Device", "factor(Time)", "Device:factor(Time)")

qBH <- matrix(p.adjust(pvals, method = "BH"), nrow = nrow(pvals))
rownames(qBH) <- colnames(df)[14:42]
colnames(qBH) <- colnames(pvals)

sigvars <- apply(apply(qBH, 2, function(x) x<=0.05), 1, function(x) sum(x) > 0)
sig.qtable <- qBH[sigvars,][order(apply(qBH[sigvars,], 1, min)),]

qtable <- as.data.frame(signif(sig.qtable, 2))
qtable %>%
    mutate(
        `B-cell` = row.names(.),
        Device = cell_spec(Device, color = ifelse(qtable$Device > 0.05, "grey", "red")),
        `Device:factor(Time)` = cell_spec(`Device:factor(Time)`, color = ifelse(qtable$`Device:factor(Time)` > 0.05, "grey", "red")),
        `factor(Time)` = cell_spec(`factor(Time)`, color = ifelse(qtable$`factor(Time)` > 0.05, "grey", "red"))
    ) %>%
    kable( escape = F,
           digits = 20,
           row.names = T,
           caption = "Significant device-type q-values") %>%
    kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
                  font_size = 12) %>%
    scroll_box(width = "100%")
Significant device-type q-values
Device factor(Time) Device:factor(Time) B-cell
1 0.16 0.071 3e-10 CD19+CD5+
2 0.33 0.041 8.7e-06 CD27+IgD- switched memory
3 0.54 0.071 7.8e-05 CD19+CD27+
4 0.54 0.31 7.8e-05 CD27-IgD- switched memory
5 0.59 4e-04 0.0017 CD19+CD11b+
6 0.071 0.16 0.0034 CD27-IgD+ mature naive
7 0.93 0.0058 0.98 CD19 of live lymph
8 0.93 0.94 0.012 CD19CD24hiCD38-memory
9 0.27 0.041 1 CD27+IgD-IgM+ switched memory

We plotted the average across time for each of the B-cells that showed a statistically significant effect across devices in the above mixed effect models. We drew attention to specific features that induced the positive test result, by listing the model parameters with \(p<0.05\) in the multivariate fit. Note that the reference level for the time comparisons is timepoint 0, and the reference level for the devices is the HeartMate-II (HMII).

require(reshape2, quietly = T)
df.long <- melt(df, id.vars = colnames(df)[1:13])
groups <- make.names(c("Device Type"))
names(df.long) <- make.names(names(df.long))

invisible(suppressMessages(require(Hmisc, quietly = T)))
stat_sum_df <- function(fun, geom="errorbar", ...) {
    stat_summary(fun.data = fun, geom = geom, width = 1, ...)
}

plots.ts <- mclapply(rownames(qtable), function(this_var){
    lapply(groups, function(this_groups){
        ggplot(subset(df.long, df.long$variable == this_var)) +
            aes(x = Time, y = value, group = PatientID) +
            aes_string(color = this_groups, fill = this_groups) +
            geom_line(alpha = 0) + 
            geom_point(alpha = 0) + 
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("point"), position = position_dodge(.5)) + 
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("line"), position = position_dodge(.5)) + 
            stat_sum_df(function(x) mean_cl_normal(x, conf.int = 0.68), mapping = aes_string(group = this_groups), position = position_dodge(.5)) + 
            #stat_smooth(aes_string(group = this_groups), method = "loess", span = 1) +
            scale_color_d3() + scale_fill_d3() +
            xlab("Time (days after surgery)") +
            ylab(this_var) + 
            ggtitle(paste(this_var)) +
            theme_classic()
    })
}, mc.cores = detectCores()-1)
names(plots.ts) <- rownames(qtable)

for(ii in c(1:length(plots.ts))){
    for(jj in 1:length(plots.ts[[ii]])){
        sumtable <- suppressMessages(summary(models.device[[rownames(qtable)[ii]]]$model))
        sumtable <- sumtable$coefficients[-1, c(4:5)] # drop intercept
        sigsum <- sumtable[sumtable[,2] <= 0.05, , drop = F]
        cat("  \n###", rownames(qtable)[ii], "\n")
        print(kable(sigsum[order(sigsum[,2]),], row.names = T) %>%
                  kable_styling(bootstrap_options = c("striped", 
                                                      "hover", 
                                                      "condensed",
                                                      "responsive"),
                                font_size = 12)
        )
        cat("  \n")
        suppressWarnings(print(plots.ts[[ii]][[jj]]))
        cat("  \n")
    }
}

CD19+CD5+

t value Pr(>|t|)
Device.TypeTAH:factor(Time)5 7.156565 0.0000000
Device.TypePVAD 3.590703 0.0008326
Device.TypePVAD:factor(Time)3 -2.877181 0.0052873
Device.TypePVAD:factor(Time)5 -2.809579 0.0063611
Device.TypePVAD:factor(Time)8 -2.720946 0.0082277

CD27+IgD- switched memory

t value Pr(>|t|)
Device.TypePVAD:factor(Time)3 -5.406064 0.0000008
Device.TypePVAD 4.233715 0.0001521
Device.TypePVAD:factor(Time)8 -3.571943 0.0006518
Device.TypePVAD:factor(Time)5 -3.404439 0.0010914
Device.TypeHVAD:factor(Time)1 2.457443 0.0165160
Device.TypeTAH:factor(Time)1 2.120071 0.0376068

CD19+CD27+

t value Pr(>|t|)
Device.TypePVAD:factor(Time)3 -4.796905 0.0000087
Device.TypePVAD 3.502922 0.0011749
Device.TypePVAD:factor(Time)5 -3.024465 0.0034502
Device.TypePVAD:factor(Time)8 -2.789363 0.0068193
Device.TypeHVAD:factor(Time)1 2.452316 0.0167262

CD27-IgD- switched memory

t value Pr(>|t|)
Device.TypeTAH:factor(Time)5 5.189516 0.0000021
Device.TypeTAH:factor(Time)14 2.262837 0.0264767
Device.TypePVAD:factor(Time)3 -2.144030 0.0352412

CD19+CD11b+

t value Pr(>|t|)
Device.TypeTAH:factor(Time)14 5.250946 0.0000015
Device.TypeTAH:factor(Time)1 2.722190 0.0081768
Device.TypePVAD:factor(Time)14 2.473834 0.0157228
Device.TypeTAH:factor(Time)3 2.083219 0.0408780

CD27-IgD+ mature naive

t value Pr(>|t|)
Device.TypeTAH:factor(Time)5 -3.905748 0.0002179
Device.TypePVAD -3.456469 0.0012664
Device.TypeTAH:factor(Time)14 -2.953963 0.0042493
Device.TypeTAH:factor(Time)1 -2.579593 0.0120122
factor(Time)8 -2.477197 0.0157031
Device.TypeHVAD:factor(Time)1 -2.083232 0.0409155

CD19 of live lymph

t value Pr(>|t|)
factor(Time)1 3.054598 0.0031913
Device.TypeHVAD:factor(Time)5 2.297061 0.0246707
Device.TypeHVAD:factor(Time)3 2.039507 0.0452529

CD19CD24hiCD38-memory

t value Pr(>|t|)
Device.TypeTAH:factor(Time)5 -3.776490 0.0003363
Device.TypeHVAD:factor(Time)1 2.344675 0.0219407
factor(Time)1 -2.161244 0.0341830

CD27+IgD-IgM+ switched memory

t value Pr(>|t|)
factor(Time)5 3.237317 0.0018346
factor(Time)3 2.819578 0.0062189
factor(Time)14 2.232742 0.0287354

HeartMate-II analysis

The HeartMate-II (HMII) recipients were the largest group, and we analyzed them by themselves due to the previously observed variability across devices.

require(reshape2, quietly = T)
df.HMII <- subset(df, df$`Device Type`=="HMII")
df.long <- melt(df.HMII, id.vars = colnames(df)[1:13])

groups <- make.names(c("AgeGreater60", 
                       "Sex",
                       "LowIntermacs",
                       "RVAD", 
                       "Sensitized",
                       "VAD Indication", 
                       #"Device Type", 
                       "Survival",
                       "Outcome"))

names(df.long) <- make.names(names(df.long))

plots.ts <- mclapply(unique(df.long$variable), function(this_var){
    lapply(groups, function(this_groups){
        ggplot(subset(df.long, df.long$variable == this_var)) +
            aes(x = Time, y = value, group = PatientID) +
            aes_string(color = this_groups, fill = this_groups) +
            geom_line(alpha = 0) + 
            geom_point(alpha = 0) + 
            #stat_smooth(aes_string(group = this_groups), method = "loess", span = 1) +
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("point"), position = position_dodge(.5)) + 
            stat_summary(fun.y = mean, aes_string(group = this_groups), geom=c("line"), position = position_dodge(.5)) + 
            stat_sum_df(function(x) mean_cl_normal(x, conf.int = 0.68), mapping = aes_string(group = this_groups), position = position_dodge(.5)) + 
            scale_color_aaas() + scale_fill_aaas() +
            xlab("Time (days after surgery)") +
            ylab(this_var) + 
            ggtitle(paste(this_var)) +
            theme_classic()
    })
}, mc.cores = detectCores()-1)
names(plots.ts) <- unique(df.long$variable)

# for(ii in c(1:length(plots.ts))){
#     for(jj in 1:length(plots.ts[[ii]])){
#         suppressWarnings(print(plots.ts[[ii]][[jj]]))
#     }
# } 

PCA

Using PCA, we found large variability between individual patients, compared to the variability within individual patients. None of the other features were clearly separable.

suppressMessages(require(ggbiplot, quietly = T))
require(ggsci, quietly = T)

# Efron's double standardization
double_standardize <- function(x, niter = 100) {
    for(i in 1:niter) x <- t(scale(t(scale(x))))
    return(as.data.frame(x))
}

isna <- unique(unlist(apply(df.HMII[,14:42], 2, function(x) which(is.na(x)))))
pca <- prcomp(double_standardize(df.HMII[-isna, 14:42]), center = TRUE, scale. = TRUE)

colorfun <- function(grouping, ...){
    if(nlevels(factor(grouping)) > 10) scale_color_discrete(...)
    if(nlevels(factor(grouping)) <=10) scale_color_d3(...)
} 

plots.pca <- mclapply(names(df.HMII)[1:17], function(this_var){
    this_groups <- df.HMII[-isna, this_var]
    if(is.numeric(this_groups)) this_groups <- cut(this_groups, 4)
    
    ggbiplot(pca,
             groups = this_groups,
             ellipse = TRUE,
             alpha = 0.3,
             varname.size = 1.2) + 
        colorfun(this_groups, name = this_var) +
        ggtitle(paste0(this_var, " variability")) +
        theme_classic() 
}, mc.cores = detectCores()-1)

names(plots.pca) <- names(df.HMII)[1:17]
#plots.pca$PatientID
#plots.pca$`Device Type`

for(ii in 1:length(plots.pca)){
    cat("  \n###", names(plots.pca)[ii], "\n")
    suppressWarnings(print(plots.pca[[ii]]))
    cat("  \n")
}

PatientID

Age

AgeGreater60

Sex

LowIntermacs

InterMACS

RVAD

Sensitized

VAD Indication

Device Type

Survival

Outcome

Time

num Total PBMC

num lymph

lymph

live lymph

Two-way repeated measures anova

We attempted to analyze the differences in B-cell levels for various features using a mixed effect model (a.k.a. two-way repeated measures ANOVA). Here we report variables that had a statistically significant variance (Benjamini-Hochberg \(q<0.05\)) across groups, or groups at each timepoint.

suppressMessages(require(lmerTest, quietly = TRUE))
suppressMessages(require(car, quietly = TRUE))
require(reshape2, quietly = TRUE)

df.lmer <- df.HMII
names(df.lmer) <- make.names(names(df.lmer), unique = TRUE)

groupvars.ix <- c(3,4,5,7,8,9,11)
groupvars <- names(df.lmer)[groupvars.ix]

bcells.ix <- c(14:42)
bcells <- names(df.lmer)[bcells.ix]

models.b <- mclapply(groupvars, function(this_groupvar){
    models.bcells <- lapply(bcells, function(this_bcell){
        this_formula <- as.formula(paste0(this_bcell, " ~ ", this_groupvar, 
                                          " * factor(Time) + (1|PatientID)"))
        suppressMessages(suppressWarnings(this_model <- lmer(this_formula, data = droplevels(df.lmer))))
        this_anova <- Anova(this_model, type = 2)
        this_pvalues <- this_anova$`Pr(>Chisq)`
        names(this_pvalues) <- rownames(this_anova)
        #return(this_pvalues)
        return(list(model = this_model,
                    pvals = this_pvalues))
    })
    names(models.bcells) <- colnames(df)[14:42]
    pvalues <- do.call(rbind, lapply(models.bcells, function(x) x$pvals))
    rownames(pvalues) <- bcells
    #return(pvalues)
    return(list(model = models.bcells,
                pvals = pvalues))
}, mc.cores = detectCores()-1)
names(models.b) <- groupvars
pvals <- lapply(models.b, function(x) x$pvals)
# something wrong here
names(pvals) <- groupvars
pvals.matrix <- do.call(cbind, lapply(pvals, function(this_pval) this_pval[,c(1,3)]))


# Benjamini Hochberg
# qBH <- matrix(p.adjust(as.numeric(pvals.matrix), 
#                        method = "BH"), 
#               nrow = nrow(pvals.matrix), 
#               ncol = ncol(pvals.matrix), 
#               byrow = F) 
# rownames(qBH) <- rownames(pvals.matrix)
# colnames(qBH) <- colnames(pvals.matrix)
# rownames(qBH) <- names(df)[bcells.ix]
# qvalsBH.df <- melt(qBH)
# colnames(qvalsBH.df) <- c("B-cell", "parameter", "qvalue")
# qvalsBH.df.ranked <- qvalsBH.df[order(qvalsBH.df$qvalue, decreasing = F),]
# qvalsBH.df.ranked[qvalsBH.df.ranked$qvalue <= 0.3,]

# Local FDR
require(fdrtool, quietly = T)
invisible(suppressMessages(fdrobj <- fdrtool(as.numeric(pvals.matrix), statistic = "pvalue", plot = F, verbose = F)))
qvals.matrix <- matrix(fdrobj$q, nrow = nrow(pvals.matrix), ncol = ncol(pvals.matrix), byrow = F)
rownames(qvals.matrix) <- rownames(pvals.matrix)
colnames(qvals.matrix) <- colnames(pvals.matrix)
rownames(qvals.matrix) <- names(df)[bcells.ix]
qvals.df <- melt(qvals.matrix)
colnames(qvals.df) <- c("B-cell", "parameter", "qvalue")
qvals.df.ranked <- qvals.df[order(qvals.df$qvalue, decreasing = F),]
shortlist <- qvals.df.ranked[qvals.df.ranked$qvalue <= 0.3,]


kable(shortlist, 
      digits = 3,
      row.names = T,
      caption = "Significant results") %>%
    kable_styling(bootstrap_options = c("striped", "hover", "condensed", "responsive"),
                  font_size = 10) %>%
    scroll_box(width = "100%")
Significant results
B-cell parameter qvalue
2 num lymph AgeGreater60 0.040
145 CD19+27+IgD-38++IgG ASC LowIntermacs 0.040
52 CD268 of +27-38++transitional AgeGreater60:factor(Time) 0.073
276 CD27+IgD-IgM+ switched memory Sensitized:factor(Time) 0.081
61 lymph Sex 0.200
348 CD19+27+IgD-38++IgG ASC VAD.Indication:factor(Time) 0.213
232 CD19+27+IgD-38++IgG ASC RVAD:factor(Time) 0.230
155 CD27-38++ transitional LowIntermacs:factor(Time) 0.254
3 lymph AgeGreater60 0.257

We plotted the average across time for each of the B-cells that showed a statistically significant effect across various factors in the above mixed effect models. We drew attention to specific features that induced the positive test result, by listing the model parameters with \(p<0.05\) in the multivariate fit.

require(stringr, quietly = T)
siggroups <- sapply(str_split(shortlist$parameter, ":"), function(x) x[1])
for(ii in 1:nrow(shortlist)){
    this_group <-siggroups[ii]
    this_bcell <- as.character(shortlist$`B-cell`[ii])
    cat("  \n###", as.character(shortlist$`B-cell`[ii]), "\n")
    
    sumtable <- suppressMessages(summary(models.b[[this_group]]$model[[this_bcell]]$model))
    sumtable <- sumtable$coefficients[-1, c(1:5)] # drop intercept
    sigsum <- sumtable[sumtable[,5] <= 0.05, , drop = F]
    
    print(kable(sigsum[order(sigsum[,5]),,drop=F], row.names = T) %>%
              kable_styling(bootstrap_options = c("striped", 
                                                  "hover", 
                                                  "condensed",
                                                  "responsive"),
                            font_size = 12)
    )
    cat("  \n")
    
    suppressWarnings(print(plots.ts[[shortlist$`B-cell`[ii]]][[which(groups == this_group)]]))
    cat("  \n")
}

num lymph

Estimate Std. Error df t value Pr(>|t|)
factor(Time)21 169798.9 44963.12 1780.755 3.776403 0.0001643
AgeGreater60older:factor(Time)21 -214015.6 62340.16 1651.364 -3.433029 0.0006116
AgeGreater60older:factor(Time)5 -136885.7 61168.57 2679.738 -2.237844 0.0253131
factor(Time)14 128214.1 57636.41 1470.433 2.224532 0.0262639

CD19+27+IgD-38++IgG ASC

Estimate Std. Error df t value Pr(>|t|)
LowIntermacsHigh 3.095 1.279721 69 2.418497 0.0182277

CD268 of +27-38++transitional

Estimate Std. Error df t value Pr(>|t|)
AgeGreater60older:factor(Time)8 -27.22177 10.33554 54.58048 -2.633803 0.0109633
factor(Time)14 -24.92941 12.43679 55.44507 -2.004489 0.0499111

CD27+IgD-IgM+ switched memory

Estimate Std. Error df t value Pr(>|t|)
SensitizedYes:factor(Time)3 17.965877 5.302290 23.31974 3.388324 0.0024949
factor(Time)5 9.747852 4.547781 23.30369 2.143431 0.0427391

lymph

Estimate Std. Error df t value Pr(>|t|)
SexMale -22.22807 9.001586 53.16972 -2.46935 0.0167811

CD19+27+IgD-38++IgG ASC

Estimate Std. Error df t value Pr(>|t|)
factor(Time)5 4.722429 1.330913 60.14685 3.548263 0.0007592
VAD.IndicationDT:factor(Time)5 -6.577445 1.874482 57.09611 -3.508940 0.0008850

CD19+27+IgD-38++IgG ASC

Estimate Std. Error df t value Pr(>|t|)
RVADYes:factor(Time)5 11.55470 2.912778 59.01744 3.966900 0.0001997
factor(Time)14 -2.17979 1.043488 56.60858 -2.088945 0.0412202

CD27-38++ transitional

Estimate Std. Error df t value Pr(>|t|)
LowIntermacsHigh:factor(Time)21 -11.380113 2.894440 59.65377 -3.931715 0.0002222
factor(Time)21 9.440964 2.450218 61.33982 3.853112 0.0002815

lymph

Estimate Std. Error df t value Pr(>|t|)